1. Load Data This markdown file contains the scRNA-seq analysis for the paper titled Integrated protein and transcriptome high-throughput spatial profiling. The data used for this analysis can be found under the GEO repository GSE198353.
# install Seurat v4.0.0 and spots
if (!requireNamespace("remotes", quietly = TRUE))
  install.packages("remotes")
if (!requireNamespace("Seurat", quietly = TRUE) | utils::packageVersion("Seurat") < "4.0.0")
  remotes::install_version("Seurat", version = "4.0.0")
if (!requireNamespace("spots", quietly = TRUE))
  install.packages("spots")

# load data
library(Seurat)
library(spots)

mmtv_gex <- Read10X_h5('GSE198353_mmtv_pymt_GEX_filtered_feature_bc_matrix.h5')
mmtv_adt <- read.csv('GSE198353_mmtv_pymt_ADT.csv.gz', header = TRUE, row.names = 1, check.names = FALSE)
mmtv_image <- Read10X_Image('GSE198353_mmtv_pymt_spatial')
mmtv <- CreateSeuratObject(mmtv_gex, assay = "RNA", project = "MMTV")
mmtv_adt <- CreateSeuratObject(mmtv_adt, assay = "CITE", project = "MMTV")
mmtv@assays$CITE <- mmtv_adt@assays$CITE
mmtv$nCount_CITE <- mmtv_adt$nCount_CITE
mmtv$nFeature_CITE <- mmtv_adt$nFeature_CITE
mmtv_image@key <- "A"
mmtv@images <- list(A = mmtv_image)
SpatialDimPlot(mmtv)

  1. Data normalization
mmtv <- NormalizeData(mmtv, assay = "RNA", verbose = FALSE)
mmtv <- NormalizeData(mmtv, assay = "CITE", verbose = FALSE)
DefaultAssay(mmtv) <- "CITE"
mmtv <- ScaleData(mmtv, verbose = FALSE)
SpatialFeaturePlot(mmtv,features = c("CD326","Podoplanin","I-A/I-E",
                                      "F4/80", "CD11c","CD4"), ncol = 3, min.cutoff = "q25")

  1. Spatial Component Analysis and Clustering
Visium.hnn.dist <- LoadData("~/Downloads", "Visium.HNN")
mmtv.hnn.dist <- VisiumHnn("~/Downloads/", Cells(mmtv))
mmtv.hnn <- HnnNeighbor(mmtv.hnn.dist, k = 37, include.self = FALSE)
mmtv.hnn.weight <- HnnWeight(mmtv.hnn$dist.mat, dist.k = 3, sigma = 1)
mmtv.sca <- SCA(X = Matrix::t(mmtv@assays$CITE@data), 
                W = mmtv.hnn.weight, 
                scaled.data = t(mmtv@assays$CITE@scale.data),
                n.eigen = 30)
mmtv@reductions[["sca"]] <- CreateDimReducObject(embeddings = mmtv.sca$X, 
                                                 loadings = mmtv.sca$rotation, 
                                                 stdev = mmtv.sca$eigenvalues, 
                                                 key = "SC_", assay = "CITE")
mmtv <- FindNeighbors(mmtv, reduction = "sca", dims = 1:15)
mmtv <- FindClusters(mmtv, resolution = 0.3)
## Modularity Optimizer version 1.3.0 by Ludo Waltman and Nees Jan van Eck
## 
## Number of nodes: 1978
## Number of edges: 72438
## 
## Running Louvain algorithm...
## Maximum modularity in 10 random starts: 0.8278
## Number of communities: 5
## Elapsed time: 0 seconds
mmtv <- RenameIdents(mmtv, '0' = 'Fibroblast-high',
                     '1' = 'Fibroblast-low',
                     '2' = 'M2-MF',
                     '3' = 'Immune-enriched (Peritumor)',
                     '4' = 'M1-MF')
mmtv@active.ident = factor(mmtv@active.ident, levels = c('Fibroblast-high', 'Fibroblast-low','M1-MF',
                                                         'M2-MF','Immune-enriched (Peritumor)'))
SpatialDimPlot(mmtv, label = TRUE, repel = TRUE)

  1. Differentially expressed ADTs
adt.markers <- FindAllMarkers(mmtv, only.pos = TRUE, logfc.threshold = 0.2, verbose = FALSE)
DoHeatmap(mmtv, adt.markers$gene, assay = "CITE", angle = 45)

  1. Differentially expressed mRNAs
DefaultAssay(mmtv) = "RNA"
mrna.markers <- FindAllMarkers(mmtv, only.pos = TRUE, verbose = FALSE)
mmtv.avg <- AverageExpression(mmtv, assays = "RNA", return.seurat = TRUE)
DoHeatmap(mmtv.avg, mrna.markers$gene, draw.lines = FALSE, assay = "RNA", angle = 45)

  1. Session info
# Session info
print(sessionInfo())
## R version 4.0.3 (2020-10-10)
## Platform: x86_64-apple-darwin17.0 (64-bit)
## Running under: macOS Big Sur 10.16
## 
## Matrix products: default
## BLAS:   /Library/Frameworks/R.framework/Versions/4.0/Resources/lib/libRblas.dylib
## LAPACK: /Library/Frameworks/R.framework/Versions/4.0/Resources/lib/libRlapack.dylib
## 
## locale:
## [1] en_US.UTF-8/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8
## 
## attached base packages:
## [1] stats     graphics  grDevices utils     datasets  methods   base     
## 
## other attached packages:
## [1] spots_0.1.0        SeuratObject_4.0.0 Seurat_4.0.1      
## 
## loaded via a namespace (and not attached):
##   [1] Rtsne_0.15            colorspace_2.0-0      deldir_0.2-3          ellipsis_0.3.2        ggridges_0.5.3        rstudioapi_0.13      
##   [7] spatstat.data_1.7-0   farver_2.0.3          leiden_0.3.6          listenv_0.8.0         remotes_2.4.2         bit64_4.0.5          
##  [13] ggrepel_0.9.0         RSpectra_0.16-0       codetools_0.2-18      splines_4.0.3         knitr_1.38            polyclip_1.10-0      
##  [19] jsonlite_1.7.2        ica_1.0-2             cluster_2.1.0         png_0.1-7             uwot_0.1.10           shiny_1.5.0          
##  [25] sctransform_0.3.2     spatstat.sparse_2.0-0 compiler_4.0.3        httr_1.4.2            Matrix_1.3-2          fastmap_1.1.0        
##  [31] lazyeval_0.2.2        limma_3.46.0          cli_3.2.0             later_1.1.0.1         htmltools_0.5.2       tools_4.0.3          
##  [37] igraph_1.2.6          gtable_0.3.0          glue_1.6.2            RANN_2.6.1            reshape2_1.4.4        dplyr_1.0.2          
##  [43] Rcpp_1.0.7            scattermore_0.7       jquerylib_0.1.3       vctrs_0.4.1           nlme_3.1-151          lmtest_0.9-38        
##  [49] xfun_0.30             stringr_1.4.0         globals_0.14.0        mime_0.9              miniUI_0.1.1.1        lifecycle_1.0.1      
##  [55] irlba_2.3.3           goftest_1.2-2         future_1.21.0         MASS_7.3-53           zoo_1.8-8             scales_1.1.1         
##  [61] spatstat.core_1.65-5  promises_1.1.1        spatstat.utils_2.1-0  parallel_4.0.3        RColorBrewer_1.1-2    yaml_2.2.1           
##  [67] reticulate_1.18       pbapply_1.4-3         gridExtra_2.3         ggplot2_3.3.3         sass_0.4.0            rpart_4.1-15         
##  [73] stringi_1.5.3         highr_0.8             rlang_1.0.2           pkgconfig_2.0.3       matrixStats_0.57.0    evaluate_0.15        
##  [79] lattice_0.20-41       ROCR_1.0-11           purrr_0.3.4           tensor_1.5            labeling_0.4.2        patchwork_1.1.1      
##  [85] htmlwidgets_1.5.3     bit_4.0.4             cowplot_1.1.1         tidyselect_1.1.0      parallelly_1.23.0     RcppAnnoy_0.0.18     
##  [91] plyr_1.8.6            magrittr_2.0.1        R6_2.5.0              generics_0.1.0        DBI_1.1.0             mgcv_1.8-33          
##  [97] pillar_1.4.7          fitdistrplus_1.1-3    survival_3.2-7        abind_1.4-5           tibble_3.0.4          future.apply_1.7.0   
## [103] hdf5r_1.3.5           crayon_1.3.4          KernSmooth_2.23-18    spatstat.geom_1.65-5  plotly_4.9.2.2        rmarkdown_2.13       
## [109] grid_4.0.3            data.table_1.13.6     digest_0.6.27         xtable_1.8-4          tidyr_1.1.2           httpuv_1.5.4         
## [115] munsell_0.5.0         viridisLite_0.3.0     bslib_0.3.1

  1. Tri-Institutional Training Program in Computational Biology and Medicine, xin2001@med.cornell.edu↩